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Abstract 

We consider equilibrium relaxation properties of the end-to-end distance and of principal com- 
ponents in a one-dimensional polymer chain model with nonlinear interaction between the beads. 
While for the single-well potentials these properties are similar to the ones of a Rouse chain, for 
the double-well interaction potentials, modeling internal friction, they differ vastly from the ones 
of the harmonic chain at intermediate times and intermediate temperatures. This minimal de- 
scription within a one-dimensional model mimics the relaxation properties found in much more 
complex polymer systems. Thus, the relaxation time of the end-to-end distance may grow by or- 
ders of magnitude at intermediate temperatures. The principal components (whose directions are 
shown to coincide with the normal modes of the harmonic chain, whatever interaction potential is 
assumed) not only display larger relaxation times but also subdiffusive scaling. 

PACS numbers: 87.15.H-, 05.20.-y, 05.40.-a 
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I. INTRODUCTION 



The dynamics of polymers and peptides attracted large attention in the past decade, 
mostly due to the relations of such dynamics to biological function. First such interest 
was due to the dynamics of nonequilibrium states (connected with the folding problem and 
with the biological functioning of the proteins), and only later some interest to equilibrium 
fluctuations in proteins arose. This interest was caused by the two reasons: one has to 
do with the thermal stability of proteins as connected to their structures [l)], another one 
is related to the luminescent measurements of fluctuations of distance between the two 
groups in equilibrium and the discovery of the anomalous kinetics and of extremely 
large characteristics times in such fluctuations. The results of these investigations, put 
together, lead to an enigma: On one hand, the thermal stability and many other properties 
of proteins can be well-described within simple bead-spring models of such systems 0, 4|, 
which for small deviations from equilibrium can be reduced to a standard picture of normal 
modes in a complex harmonic network. The evaluation of the characteristic times in such 
systems, however, leads to results being off by orders of magnitudes when compared to the 
observed times p]. The strong discrepancy in correlation times implies the existence of 
a strong additional internal friction mechanism slowing down the dynamics compared the 
linear Rouse one. 

Recently, the anomalous kinetics (power-law decay of correlations) was observed even in 
single modes and even in relatively short and flexible peptides, i.e. linear chains lacking 



secondary structure [6| , and the behavior found here is very close to the one observed in 
protein simulations [J. Therefore such kinetics probably is an intrinsic property of a large 
class of polymers, and is not connected to the specific properties of the secondary structure of 
proteins. Furthermore it has been shown that trap models cannot account for the anomalous 
dynamics [6i]. Therefore, it is necessary to analyze, what kind of the minimal assumptions 
about the intramolecular potentials have to be made to build a model mimicking the behavior 
observed in realistic molecular dynamics simulations. 

Typical mechanisms of internal friction involve the existence of barriers in the free energy 
landscape of the system, being of entropic or of enthalpic nature. The necessity to overcome 
such barriers slows down the overall dynamics by the Arrhenius factor which might be quite 
large [8;]. However, the existence of a complex, nonlinear energy landscape may lead to a 
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strong mixing of modes appearing in the linearized description, and make the whole analysis 
based on such a picture problematic. As we proceed to show, this is not the case: although 
the dynamics of modes is strongly influenced, the directions of the principal components 
(PCs) in configuration space follow those of the normal modes, which resolves the enigma. 

In what follows we first consider the Brownian dynamics of a three-dimensional polymer 
chain, as applied e.g. for simulating polyethylene molecule within the valent-angle model 
[sl, [lo|, III, and show that it follows slow kinetics at intermediate times, provided that the 
temperature is low enough. On the other hand, at high temperatures the typical Rouse 
dynamics sets on. 

To simplify the model even more, we consider a one-dimensional chain of beads, similar 
to the Rouse model of the polymer dynamics, with the only difference being changing the 
harmonic interaction between the beads for a nonlinear one. We discuss the end-to-end 
distance of the chain and the relaxation properties of the PCs 12| of the system. 

The one- dimensional model discussed here is a close continuous analogue of the so-called 
"necklace model" for reptation 13|, llJ] , so that the results might have a broader applicability 
also outside of the scope of present investigations. 

The results of extensive numerical simulations show that while a whatever single-well po- 
tential tested does not affect drastically single-mode relaxation properties, the double-well 
potentials lead to intermediate-time behavior quite similar to the one observed in realistic 
systems and to strong increase in relaxation times compared to the harmonic case. Therefore 
such a model may be considered as a possible candidate for the explanation of the corre- 
sponding findings. The multiwell intramolecular potentials appear at different scales and 
are ubiquitous in polymers [sl, Q. They appear quite naturally e.g. within the valent-angle 
model. The anomalous behavior is present within the finite temperature range (which in the 
realistic case may be relevant for the biological functioning), and disappears both for low 
and for high temperatures, where the dynamics can therefore be described by an effective 
harmonic model, albeit with the values of parameters strongly different from the "micro- 
scopic" ones. This property can explain that the generalized Gaussian models work quite 
well in predicting thermal stability properties while failing to mimic the temporal fluctuation 
behavior at moderate temperatures. 

The work is organized as follows: In the next section we present the relaxation properties 
found in a complex three-dimensional polymer model and introduce in Sec. IIIII the one- 
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dimensional chain model. In Sec. [IV]we consider the relaxation properties of the end-to-end 
distance of the chain and discuss in Sec. |V] the behavior of its PCs. We focus on their 
directions and their temporal scaling, which is found to coincide with the one found in the 
realistic three-dimensional model. Finally we summarize our results. 



II. INTERNAL FRICTION IN A THREE-DIMENSIONAL CHAIN MODEL 



In the first part of the work we consider the relaxation properties of a complex polymer 
model which was proposed to mimic a polyethylene molecule. The model equations can 



be found 



Lien was p: 
in 0, Q 



They include valence bond-, valence angle- and torsional angle- 
interactions. For our purpose we neglect the Lenard Jones interactions and focus on the 
influence of the angle interactions which follow a multistable potential energy landscape. 
All parameter values are the same as in [11] except for the bond length Iq = 1., and the 
constants ki, = 3., = 3., and /c^ = 0.1 which only set the timescale of simulation. We 
perform Brownian dynamics simulations and study the relaxation properties of the PCs and 
their autocorrelation functions. 

For a discrete set of M measurements of a trajectory x{t) the unbiased autocorrelation 
function (ACF) of an observable x{t) (coordinate, end-to-end distance, etc.) 0(t) (t G 
0, M — 1) is given by 

^ M-t-l 

4>{t) = , {<n + t)-^l){x{n)-^l) , (1) 

vvi i}a 

with the sample mean yU and sample variance cr^. Since (for /i = 0) the mean squared 
displacement {x'^it)) = {{x{t + 1') - x{t')f) = {x^{t + t') - 2x(t + t')x{t') + x^{t')) = 
2cr^ (1 — 0(t)), the values of 1 — (pit) and {x'^if)) (e.g. used in [Q]) contain the same in- 
formation and differ only in normalization. We denote the unbiased ACF of the unbiased 
ACF of the fc-th PC by 0^. 

In Fig. [1] we show the scaling of 1 — 0jt(t) for some of the PCs. In panel (a) we have chosen 
a high value of the temperature yielding a trans-gauche barrier height of AEtg/kBT = 0.66, 
a cis barrier height of AEc/kBT = 2.4 and AEt^g/kBT = 0.16. All PCs show a scaling 

with a = 1, which is the same as for the relaxation of single modes in the Rouse chain 
Q. The scaling fo. a .owe. value of the temperature w.th AE.Jt.T ^ 7.6. a c,s barner 
height of AEc/kBT = 27.6 and AEt^g/kBT = 1.8 is presented in panel (b). The relaxation 



behavior strongly differs compared to the previous case. For small times 1 — (pkit) still 
follows follows the scaling (dotted line); but at longer time it crosses over to anomalous, 
subdiffusive behavior. The curve's slope is slightly higher for the first mode. Compared to 
the high temperature case, the typical relaxation times are shifted to longer times. Thus 
the existence of barriers in the free energy landscape at intermediate temperatures slows 
down the dynamics. However, the theoretical analysis of a three-dimensional model is still 
too involved, so that we simplify it even more in the next section. 



III. THE MINIMAL MODEL 



In order to mimic the behavior found in the complex three-dimensional model we consider 
a one dimensional chain of beads with coordinates qi, qn- The interaction potential is 
given by W{li), with /j = — being the relative displacements of the neighboring beads, 
i = 1, - ■■ ,N. The chain is free, hence W{1()) = W{l]\f^i) = 0. For a harmonic potential 
W{li) this model corresponds to a Rouse chain in Id. The dynamics is overdamped and the 
system of coupled Langevin equations reads 



= {w{h^,) + W{k)} + W2^e. , (2) 

with being Gaussian white noise of strength ksT/'y. The friction coefficient 7 is set to 
unity in what follows. The simulations correspond to solving the set of coupled equations ([2]) 
by use of a Heun integration scheme. In order to get reliable statistics numerical simulations 
have to cover many orders of magnitude in time. Therefore we confine ourselves to relatively 
short chains consisting of = 25 beads. 

We concentrate on three simple prototypes of coupling potentials W{li), either single- 



coupling potentials we consider a soft and asym- 



161 ] and a hard and symmetric quartic interaction 



well or double- well. For the single- wel 
metric Toda interaction (T-potential) 
(Q-potential), i.e, W{li) = + exp(— /j) — 1 and W{li) = \lf + respectively. The double- 
well potential (DW-potential) has the form W{li) = ^If — ^If. The latter potential exhibits 
two minima separated by a maximum of height AE = 6^/(4a). We define e = AE/{kBT). 
Note, that for small deviations from their equilibria all potentials can be approximated by 
a harmonic spring with coupling constant k, = 1 for the T-and Q-potentials and k, = 2b for 
the DW-potential. 



IV. RELAXATION OF THE END-TO-END DISTANCE 



First, let us consider the ACF of the end-to-end distance Rete = In — Qi, see Fig. [2l We 
denote the unbiased ACF of the end-to-end distance by (pete- Panel (a) shows the scaling 
of 1 — (pete vs. time for the different single-well potentials. Compared to the corresponding 
harmonic chain, the autocorrelation time reduces for the Q-potential while it grows for 
the T-potential, and this effect is stronger for higher temperatures. An explanation is 
as follows: the relaxation time of Rete is to a good approximation the largest relaxation 
time of the harmonic chain. This one behaves as ri ~ 1/k [isl. For a hard Q-potential 
the effective value of k is larger than for a harmonic one and grows with temperature, so 
that Ti is expected to be smaller, while for the soft T-potential the opposite is true. In 
panel (b) 1 — (pete is depicted for the DW-potential. For different values of e the curves 
strongly differ. At a low temperature (e = 20, dotted line) the curve coincides with the 
one of the corresponding harmonic chain (solid gray line) for a finite simulation time (the 
beads hardly jump between the two minima of the coupling potential): this is exactly what 



should be observed when such a conformation change mechanism is frozen out [17|. For 
intermediate barrier heights (e = 10, dashed-dotted line) the behavior of the ACF is strongly 
different. It shows distinctly different behaviors for short and long times interpolated by 
a plateau. The characteristic relaxation time in the long-time domain (where the kinetics 
is dominated by rare fluctuations following the Arrhenius law) is by about four orders of 
magnitude larger than for the harmonic chain. We also considered a case when the barrier 
heights are randomly distributed with the density p(e) ~ exp(A(e — emin)), e > e^m = 5, 
A = 0.1 with a cut-off at Cmax = 11 to avoid freezing (solid black line). In this case the 
plateau is smoothed out, but the increase of characteristic relaxation time persists. For 
high temperatures (e = 1, dashed line) the correlation time is much smaller compared to 
the previous case (but still larger than in the harmonic case due to the softer potential) 
and its scaling follows the one of the harmonic chain. We conclude that the increase of 
the correlation time (related to the timescale of rare barrier crossings) becomes large at 
intermediate temperatures. Qualitatively the observed behavior is in full agreement with 
the one made for the realistic three-dimensional model. 
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V. PRINCIPAL COMPONENTS (PCS) 



A. Directions of the PCs 

Let us now turn to the behavior of principal components. We proceed to show that in a 
homogeneous hnear chain with symmetric interaction potentials W{li) between the beads, 
with the overall potential energy given hy W — the PCs follow the directions of 

the normal modes of the harmonic chain, independently on the exact form of the potential 
W{li). This is due to the fact that the directions of PCs are essentially thermodynamical 
quantities. 

The proof of this fact involves the following steps. Consider a grafted chain oi N + 1 
beads, whose first bead is attached to the origin of coordinates by a weak spring, with 
the small elastic constant e (an "asymptotically free" chain). Let q be the vector of the 
coordinates, and 1 be the vector with components Iq = Qq and li for i > 1. For a given 
interaction potential the distribution of 1 is given by 

p{lo,...,ln) = ^exp 

This overall canonical distribution factorizes in the product of the distributions of single 
Zj-components, which means that the corresponding random variables are independent. Due 
to symmetry of all interaction potentials, the mean values of k are zero. Hence, the corre- 
sponding variables are uncorrelated too: (/q) = Iq = 2kBT/e, (If) = for i = 1, ...,N and 
(hlj) = for all i ^ j- The variables Qi are obtained through /j via linear transformation, 
Qi — Yl]=o^jj ^^^^ their correlator Cij — {qiQj) = + niin{i,j}Z^. The matrix C can 
be diagonalized. Its eigenvectors are the PCs of the grafted chain, whether harmonic or 
not. The normalized PCs can not depend on P or 1^ independently, but only on their ratio 
f/ll = t^e/ {2kBT). Note that for the free chain, e ^ Iq ^ oo, the PCs of a harmonic 
and anharmonic chain coincide. The same is true also for asymmetric interaction potentials, 
where the lengths li have to be corrected for thermal expansion. 

Turning to a harmonic chain, we can change U for in Eq.(l), which now reads 

p{Qq,-;Qn) = ;^exp 

where M is the tridiagonal force matrix, whose elements are: moo — i + s, rriNN — 1, the 
rest of diagonal elements ma — 2, the elements on the sub- and super-diagonals are equal 
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2kBT 



qMq 



(4) 



to —1. Since the corresponding distribution is a multivariate Gaussian, the elements of the 
correlation matrix C are proportional to the ones of the inverse of the matrix M. The 
eigenvectors of the matrix M are the normal modes of the harmonic chain. Since the matrix 
M and its inverse C shear the eigenvectors, those are also the PCs of the harmonic chain. 
The idea of considering the grafted chain is connected to the wish to have M invertible. In 
the last step one shows that for e — the limits of all corresponding eigenvectors exist (the 
fact having to do with the nondegenerate nature of a harmonic chain's spectrum). 

Thus assuming whatever interaction potential between neighboring beads, the PCs follow 
the directions of the normal modes of the harmonic chain. However, as we proceed to show, 
the dynamical properties of each single mode can be drastically different from the ones of a 
harmonic chain. 

B. Relaxation properties of the PCs 

From the normal mode analysis of a harmonic chain it is known that the autocorrelation 
time of the k—th normal mode scales as ~ N'^/k'^ [isl, [3] and does not depend on the 
noise strength ksT/'-f. Although the autocorrelation time becomes temperature dependent 
for the single- well nonlinear potentials, the scaling of 1 — remains the same as in the 
harmonic situation. This is shown in Fig. [3] in panel (a). The curves order relative to the 
corresponding harmonic limiting case is the same as in Fig. O Furthermore we find that 
the scaling with k and N persits (not shown). Thus as in the harmonic chain higher modes 
have smaller correlation times. 

In contrast, for the DW-potential with intermediate barrier height of e = 10, Fig. [3] panel 
(b), the relaxation times of the lowest three modes are almost equal: the corresponding 
curves merge for t > 10^. The overall shape of the curves is similar to the one for Rete, 
albeit with different scaling. For a harmonic chain the scaling is given by with a = 1 
(corresponding to the dotted line in the plot). For small times 1 — (pkit) still follows follows 
this scaling; but at longer time it crosses over to anomalous, subdiffusive behavior with 
a ~ 0.82 over two orders of magnitude in time. The case of distributed barrier heights is 
shown in panel (c). Here we observe again a subdiffusive scaling of the ACF. Interestingly 
the shown curves for different values of k do not merge at longer times and their slope is 
mode dependent with smaller a for higher values of k {k = 1: a ~ 0.79, k = 2: a ~ 
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0.72, k = 3: a ~ 0.68). For modes k > 4 the slope changes only slightly (not shown). 
The relaxation behavior of the realistic chain presented in Fig. [1] is the same, albeit with 
somewhat different exponents. However, the exponents are expected to depend crucially on 
the choice of parameters and the ratio of different timescales in the system. The proposed 
minimal model is able to mimic qualitatively the relaxation behavior of a much more complex 
structure and therewith it helps to understand the underlying mechanism of internal friction. 

VI. SUMMARY 

To summarize, it was found that under nonlinear interaction potentials, the scaling of the 
end-to-end ACF is essentially the one of the harmonic chain, though with shifted correlation 
times. While these changes are not too large for the single- well potentials, for the double- well 
ones (describing the internal friction) these can lead to the increase in the relaxation times 
by the orders of magnitude. Turning to the principal components we have shown that for 
homogeneous chains they follow the normal modes of the harmonic chain, but can exhibit 
vastly different kinetics, both with respect to the characteristic times and to the overall 
scaling. For the double-well potentials the last can correspond to subdiffusion. It was shown 
that the same behavior can be found in realistic polymer chain models. Therefore a one- 
dimensional chain with double-well potentials might be a possible minimal model describing 
the similar observations of anomalous kinetics in realistic molecular dynamics simulations. 
Moreover, the model is a close continuous analogue of the so-called "necklace" model for 
reptation, so that the results might have a broader apphcability also outside of the scope of 
present investigations. 
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Figure captions 



Figure [T]. 

l — (f)k{t) of some of the PCs. Panel (a): High temperature value. Panel (b): Intermediate 
temperature value for which strong internal friction is observed. The parameter values are 
given in pjj] and in the text. 

Figure ^ 

ACF of the end-to-end distance Rete- Shown is 1 — (peteif)- Panel (a): Single- well inter- 
action potentials. The noise strengths is ksT/'y = 5. Panel (b): DW-potential. The barrier 
heights are given in the legend, 6 = 2. 

Figure [S]. 

(Color online) 1 — (f)k{t) of the lowest PCs. Panel (a): Single- well potentials. The noise 
strength is the same as in Fig. [2l Panel (b): DW-potential with the parameter values a = 1, 
6 = 2, and e = 10. Panel (c): DW-potential with barrier heights distributed as described in 
the text. 
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